Functions
produce_formants <- function(target, bias = 0) {
noise <- sample(seq(-0.2, 0.2, by = 0.05), 1)
outcome <- target + noise + bias
return(outcome)
}
create_lexicon <- function(lexicon_size, vowels, group_size, formants) {
lexicon <- tibble(
word = 1:lexicon_size,
vowel = rep(vowels, len = lexicon_size),
class = sample(rep(1:(lexicon_size/group_size), each = group_size,
len = lexicon_size)),
frequency = rep(1:lexicon_size, each = length(vowels), len = lexicon_size),
f1 = rep(as.list(formants), len = lexicon_size)
)
return(lexicon)
}
create_lexicon_2 <- function(lexicon_size, vowels, group_size, formants_1, formants_2) {
lexicon <- tibble(
word = 1:lexicon_size,
vowel = rep(vowels, len = lexicon_size),
class = sample(rep(1:(lexicon_size/group_size), each = group_size,
len = lexicon_size)),
frequency = rep(1:lexicon_size, each = length(vowels), len = lexicon_size),
f1 = rep(as.list(formants_1), len = lexicon_size),
f2 = rep(as.list(formants_2), len = lexicon_size)
)
return(lexicon)
}
resample <- function(x) {
if (length(x) == 1) {
return(x)
} else {
sample(x, 1)
}
}
populate_lexicon <- function(lexicon) {
lexicon_size <- nrow(lexicon)
lexicon_frequency <- lexicon$frequency
lexicon_f1 <- lexicon$f1
for (i in 1:50000) {
word_id <- sample(1:lexicon_size, 1, prob = lexicon_frequency)
target <- resample(unlist(lexicon_f1[[word_id]]))
outcome <- produce_formants(target)
lexicon_f1[[word_id]][[length(lexicon_f1[[word_id]]) + 1]] <- outcome
}
lexicon$f1 <- lexicon_f1
return(lexicon)
}
populate_lexicon_2 <- function(lexicon) {
lexicon_size <- nrow(lexicon)
lexicon_frequency <- lexicon$frequency
lexicon_f1 <- lexicon$f1
lexicon_f2 <- lexicon$f2
for (i in 1:50000) {
word_id <- sample(1:lexicon_size, 1, prob = lexicon_frequency)
target_1 <- resample(unlist(lexicon_f1[[word_id]]))
outcome_1 <- produce_formants(target_1)
lexicon_f1[[word_id]][[length(lexicon_f1[[word_id]]) + 1]] <- outcome_1
target_2 <- resample(unlist(lexicon_f2[[word_id]]))
outcome_2 <- produce_formants(target_2)
lexicon_f2[[word_id]][[length(lexicon_f2[[word_id]]) + 1]] <- outcome_2
}
lexicon$f1 <- lexicon_f1
lexicon$f2 <- lexicon_f2
return(lexicon)
}
get_encode <- function(condition) {
if (condition) {
encoding_prob <- lexicon_frequency[word_id] /
max(lexicon_frequency)
encode <- sample(c(TRUE, FALSE), 1,
prob = c(
encoding_prob, 1 - encoding_prob
)
)
} else {
encode <- TRUE
}
return(encode)
}
sound_shift <- function(lexicon_size, vowels, group_size, formants, biased_vowel, bias, iterations, save_freq) {
lexicon <- create_lexicon(lexicon_size, vowels, group_size, formants) %>%
populate_lexicon()
lexicon_size <- nrow(lexicon)
lexicon_vowel <- lexicon[["vowel"]]
lexicon_class <- lexicon[["class"]]
lexicon_frequency <- lexicon[["frequency"]]
lexicon_f1 <- lexicon[["f1"]]
new_lexicon_f1 <- tibble(init = 1:lexicon_size)
environment(get_encode) <- environment()
for (iteration in 1:iterations) {
word_id <- sample(1:lexicon_size, 1, prob = lexicon_frequency)
vowel <- lexicon_vowel[word_id]
word_class <- lexicon_class[word_id]
if (vowel == biased_vowel) {
current_bias <- bias
} else {
current_bias <- 0
}
#### Produce the chosen word ####
target <- resample(unlist(lexicon_f1[word_id]))
outcome <- produce_formants(target, current_bias)
if (vowel == vowels[1]) {
pool_max <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
encode <- get_encode(outcome <= pool_max)
} else if (vowel == vowels[2]) {
if (length(vowels) == 2) {
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)
))
encode <- get_encode(outcome >= pool_min)
} else {# if length(vowels) == 3
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)
))
pool_max <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[3])]
)
))
encode <- get_encode(outcome >= pool_min || outcome <= pool_max)
}
} else {# if vowel == vowels[3]
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)
))
encode <- get_encode(outcome >= pool_min)
}
#### Encode ####
if (encode) {
lexicon_f1[[word_id]][[length(lexicon_f1[[word_id]]) + 1]] <- outcome
}
if (iteration %% save_freq == 0) {
column <- as.character(iteration)
new_lexicon_f1 <- mutate(
new_lexicon_f1,
!!column := lexicon_f1
)
}
}
lexicon <- cbind(lexicon, new_lexicon_f1) %>%
rename(`0` = f1) %>%
select(-init) %>%
gather(time, f1, matches("\\d")) %>%
mutate(time = as.integer(time)) %>%
group_by(time, word) %>%
mutate(f1_mean = mean(unlist(f1))) %>%
ungroup()
return(lexicon)
}
sound_shift_2 <- function(lexicon_size, vowels, group_size, formants_1, formants_2, biased_vowel, bias, iterations, save_freq) {
lexicon <- create_lexicon_2(lexicon_size, vowels, group_size, formants_1, formants_2) %>%
populate_lexicon_2()
lexicon_size <- nrow(lexicon)
lexicon_vowel <- lexicon[["vowel"]]
lexicon_class <- lexicon[["class"]]
lexicon_frequency <- lexicon[["frequency"]]
lexicon_f1 <- lexicon[["f1"]]
lexicon_f2 <- lexicon[["f2"]]
new_lexicon_f1 <- tibble(init_1 = 1:lexicon_size)
new_lexicon_f2 <- tibble(init_2 = 1:lexicon_size)
environment(get_encode) <- environment()
for (iteration in 1:iterations) {
word_id <- sample(1:lexicon_size, 1, prob = lexicon_frequency)
vowel <- lexicon_vowel[word_id]
word_class <- lexicon_class[word_id]
if (vowel == biased_vowel) {
current_bias <- bias
} else {
current_bias <- 0
}
#### Produce the chosen word ####
target <- resample(unlist(lexicon_f1[word_id]))
outcome <- produce_formants(target, current_bias)
target_2 <- resample(unlist(lexicon_f2[word_id]))
outcome_2 <- produce_formants(target_2)
if (vowel == vowels[1]) {
pool_max <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
pool_max_2 <- suppressWarnings(max(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
pool_min_2 <- suppressWarnings(min(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
encode <- get_encode(
(outcome >= pool_min & outcome <= pool_max) &
(outcome_2 >= pool_min_2 & outcome_2 <= pool_max_2)
)
} else if (vowel == vowels[2]) {
if (length(vowels) == 2) {
pool_max <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_max_2 <- suppressWarnings(max(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_min_2 <- suppressWarnings(min(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
encode <- get_encode(
(outcome >= pool_min & outcome <= pool_max) &
(outcome_2 >= pool_min_2 & outcome_2 <= pool_max_2)
)
} else {# if length(vowels) == 3
pool_max <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_max_2 <- suppressWarnings(max(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_min_2 <- suppressWarnings(min(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[1])]
)))
pool_max_3 <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[3])]
)))
pool_min_3 <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[3])]
)))
pool_max_2_3 <- suppressWarnings(max(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[3])]
)))
pool_min_2_3 <- suppressWarnings(min(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[3])]
)))
encode <- get_encode(
(outcome >= pool_min & outcome <= pool_max) &
(outcome_2 >= pool_min_2 & outcome_2 <= pool_max_2) ||
(outcome >= pool_min_3 & outcome <= pool_max_3) &
(outcome_2 >= pool_min_2_3 & outcome_2 <= pool_max_2_3)
)
}
} else {# if vowel == vowels[3]
pool_max <- suppressWarnings(max(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
pool_min <- suppressWarnings(min(unlist(
lexicon_f1[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
pool_max_2 <- suppressWarnings(max(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
pool_min_2 <- suppressWarnings(min(unlist(
lexicon_f2[which(lexicon_vowel[which(lexicon_class == word_class)] == vowels[2])]
)))
encode <- get_encode(
(outcome >= pool_min & outcome <= pool_max) &
(outcome_2 >= pool_min_2 & outcome_2 <= pool_max_2)
)
}
#### Encode ####
if (encode) {
lexicon_f1[[word_id]][[length(lexicon_f1[[word_id]]) + 1]] <- outcome
lexicon_f2[[word_id]][[length(lexicon_f2[[word_id]]) + 1]] <- outcome_2
}
if (iteration %% save_freq == 0) {
column_1 <- paste("f1", as.character(iteration), sep = "_")
column_2 <- paste("f2", as.character(iteration), sep = "_")
new_lexicon_f1 <- mutate(
new_lexicon_f1,
!!column_1 := lexicon_f1
)
new_lexicon_f2 <- mutate(
new_lexicon_f2,
!!column_2 := lexicon_f2
)
}
}
lexicon <- cbind(lexicon, new_lexicon_f1, new_lexicon_f2) %>%
rename(`f1_0` = f1, `f2_0` = f2) %>%
select(-init_1, -init_2) %>%
gather(time, formant, matches("f[12]_\\d")) %>%
separate(time, c("formant_number", "time")) %>%
spread(formant_number, formant) %>%
mutate(time = as.integer(time)) %>%
group_by(time, word) %>%
mutate(f1_mean = mean(unlist(f1)), f2_mean = mean(unlist(f2))) %>%
ungroup()
return(lexicon)
}
Two vowels, F1, F2 (3)
Simulation
set.seed(seed)
lexicon_3 <- sound_shift_2(
lexicon_size = 100,
vowels = c("BART", "BAT"),
group_size = 5,
formants_1 = c(6.5, 5.5),
formants_2 = c(12.7, 13),
biased_vowel = "BART",
bias = -0.3,
iterations = 200000,
save_freq = 10000
) %>%
mutate(
freq_bin = ifelse(
frequency < max(frequency)/2,
"low",
"high"
)
)
Plotting
lexicon_3 %>%
ggplot(aes(time, f1_mean)) +
geom_jitter(alpha = 0.3, size = 0.5) +
geom_smooth(aes(colour = vowel))
## `geom_smooth()` using method = 'gam'

lexicon_3 %>%
ggplot(aes(time, f1_mean, colour = frequency, group = word)) +
geom_line() +
facet_grid(. ~ vowel)

lexicon_3 %>%
ggplot(aes(time, f1_mean, colour = freq_bin, group = word)) +
geom_line(alpha = 0.5) +
facet_grid(. ~ vowel)

lexicon_3 %>%
ggplot(aes(time, f1_mean, colour = freq_bin)) +
geom_point(alpha = 0.1, size = 0.5) +
geom_smooth(se = FALSE) +
facet_grid(. ~ vowel)
## `geom_smooth()` using method = 'loess'

lexicon_3 %>%
filter(time == 0) %>%
ggplot(aes(frequency, f1_mean)) +
geom_point(alpha = 0.1) +
geom_smooth(aes(colour = vowel), method = "lm") +
ylim(5, 7)

lexicon_3 %>%
filter(time == 100000) %>%
ggplot(aes(frequency, f1_mean)) +
geom_point(alpha = 0.1) +
geom_smooth(aes(colour = vowel), method = "lm") +
ylim(5, 7)

lexicon_3 %>%
filter(time == 200000) %>%
ggplot(aes(frequency, f1_mean)) +
geom_point(alpha = 0.1) +
geom_smooth(aes(colour = vowel), method = "lm") +
ylim(5, 7)

lexicon_3 %>%
filter(time == 0) %>%
ggplot(aes(f1_mean, f2_mean, colour = vowel)) +
geom_point() +
xlim(5, 7) + ylim(12, 13.5)

lexicon_3 %>%
filter(time == 100000) %>%
ggplot(aes(f1_mean, f2_mean, colour = vowel)) +
geom_point() +
xlim(5, 7) + ylim(12, 13.5)

lexicon_3 %>%
filter(time == 200000) %>%
ggplot(aes(f1_mean, f2_mean, colour = vowel)) +
geom_point() +
xlim(5, 7) + ylim(12, 13.5)
## Warning: Removed 1 rows containing missing values (geom_point).

Analysis
Two vowels, F1
lexicon <- lexicon %>%
mutate(
vowel_ord = ordered(vowel, levels = c("BART", "BAT")),
word_ord = as.ordered(word)
) %>%
arrange(word, time) %>%
create_event_start("word")
lexicon_gam <- bam(
f1_mean ~
frequency +
vowel_ord +
s(time, bs = "cr") +
s(time, bs = "cr", by = frequency) +
s(time, bs = "cr", by = vowel_ord) +
s(frequency, bs = "cr") +
ti(time, frequency) +
ti(time, frequency, by = vowel_ord) +
s(time, word_ord, bs = "fs"),
data = lexicon,
method = "fREML"
)
summary(lexicon_gam)
##
## Family: gaussian
## Link function: identity
##
## Formula:
## f1_mean ~ frequency + vowel_ord + s(time, bs = "cr") + s(time,
## bs = "cr", by = frequency) + s(time, bs = "cr", by = vowel_ord) +
## s(frequency, bs = "cr") + ti(time, frequency) + ti(time,
## frequency, by = vowel_ord) + s(time, word_ord, bs = "fs")
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.8451582 0.0280640 208.28 <2e-16 ***
## frequency -0.0017286 0.0009714 -1.78 0.0754 .
## vowel_ord.L -0.5890718 0.0195608 -30.11 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(time) 7.363e+00 7.524e+00 13.147 < 2e-16 ***
## s(time):frequency 5.255e+00 5.823e+00 97.520 < 2e-16 ***
## s(time):vowel_ordBAT 8.658e+00 8.823e+00 66.936 < 2e-16 ***
## s(frequency) 1.162e-04 1.165e-04 0.418 0.994
## ti(time,frequency) 7.312e+00 7.407e+00 6.534 1.84e-07 ***
## ti(time,frequency):vowel_ordBAT 7.135e+00 7.237e+00 9.155 1.96e-11 ***
## s(time,word_ord) 6.818e+02 9.930e+02 4847.119 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Rank: 1069/1072
## R-sq.(adj) = 1 Deviance explained = 100%
## fREML = -7684.5 Scale est. = 7.2608e-06 n = 2100
fvisgam(lexicon_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BART"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 50.000000.
## * vowel_ord : factor; set to the value(s): BART.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

fvisgam(lexicon_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BAT"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 50.000000.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

Three vowels, F1
lexicon_2 <- lexicon_2 %>%
mutate(
vowel_ord = ordered(vowel, levels = c("BART", "BAT", "BET")),
word_ord = as.ordered(word)
) %>%
arrange(word, time) %>%
create_event_start("word")
lexicon_2_gam <- bam(
f1_mean ~
frequency +
vowel_ord +
s(time, bs = "cr") +
s(time, bs = "cr", by = frequency) +
s(time, bs = "cr", by = vowel_ord) +
s(frequency, bs = "cr") +
ti(time, frequency) +
ti(time, frequency, by = vowel_ord) +
s(time, word_ord, bs = "fs"),
data = lexicon_2,
method = "fREML"
)
summary(lexicon_2_gam)
##
## Family: gaussian
## Link function: identity
##
## Formula:
## f1_mean ~ frequency + vowel_ord + s(time, bs = "cr") + s(time,
## bs = "cr", by = frequency) + s(time, bs = "cr", by = vowel_ord) +
## s(frequency, bs = "cr") + ti(time, frequency) + ti(time,
## frequency, by = vowel_ord) + s(time, word_ord, bs = "fs")
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.434199 0.027929 194.575 < 2e-16 ***
## frequency -0.006797 0.001431 -4.749 2.25e-06 ***
## vowel_ord.L -1.257763 0.023874 -52.683 < 2e-16 ***
## vowel_ord.Q -0.074506 0.023995 -3.105 0.00194 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(time) 7.857e+00 7.929e+00 91.547 < 2e-16 ***
## s(time):frequency 1.601e+00 1.733e+00 364.850 < 2e-16 ***
## s(time):vowel_ordBAT 8.714e+00 8.866e+00 88.733 < 2e-16 ***
## s(time):vowel_ordBET 8.675e+00 8.846e+00 78.286 < 2e-16 ***
## s(frequency) 3.792e-05 3.803e-05 0.074 0.99867
## ti(time,frequency) 1.000e+00 1.000e+00 105.938 < 2e-16 ***
## ti(time,frequency):vowel_ordBAT 1.000e+00 1.000e+00 6.761 0.00941 **
## ti(time,frequency):vowel_ordBET 1.651e+00 1.663e+00 6.107 0.01626 *
## s(time,word_ord) 6.623e+02 9.900e+02 4927.740 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Rank: 1095/1098
## R-sq.(adj) = 1 Deviance explained = 100%
## fREML = -7846 Scale est. = 6.4822e-06 n = 2100
fvisgam(lexicon_2_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BART"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 34.000000.
## * vowel_ord : factor; set to the value(s): BART.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

fvisgam(lexicon_2_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BAT"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 34.000000.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

fvisgam(lexicon_2_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BET"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 34.000000.
## * vowel_ord : factor; set to the value(s): BET.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

plot_smooth(
lexicon_2_gam,
view = "time",
cond = list(frequency = 10, vowel_ord = "BAT"),
rug = F,
col = "red",
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; set to the value(s): 10.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##
plot_smooth(
lexicon_2_gam,
view = "time",
cond = list(frequency = 20, vowel_ord = "BAT"),
rug = F,
col = "blue",
add = T,
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; set to the value(s): 20.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##
plot_smooth(
lexicon_2_gam,
view = "time",
cond = list(frequency = 30, vowel_ord = "BAT"),
rug = F,
col = "green",
add = T,
rm.ranef = TRUE
)

## Summary:
## * frequency : numeric predictor; set to the value(s): 30.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##
Two vowels, F1, F2
lexicon_3 <- lexicon_3 %>%
mutate(
vowel_ord = ordered(vowel, levels = c("BART", "BAT", "BET")),
word_ord = as.ordered(word)
) %>%
arrange(word, time) %>%
create_event_start("word")
lexicon_3_gam <- bam(
f1_mean ~
frequency +
vowel_ord +
s(time, bs = "cr") +
s(time, bs = "cr", by = frequency) +
s(time, bs = "cr", by = vowel_ord) +
s(frequency, bs = "cr") +
ti(time, frequency) +
ti(time, frequency, by = vowel_ord) +
s(time, word_ord, bs = "fs"),
data = lexicon_3,
method = "fREML"
)
summary(lexicon_3_gam)
##
## Family: gaussian
## Link function: identity
##
## Formula:
## f1_mean ~ frequency + vowel_ord + s(time, bs = "cr") + s(time,
## bs = "cr", by = frequency) + s(time, bs = "cr", by = vowel_ord) +
## s(frequency, bs = "cr") + ti(time, frequency) + ti(time,
## frequency, by = vowel_ord) + s(time, word_ord, bs = "fs")
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.8603978 0.0257960 227.183 < 2e-16 ***
## frequency -0.0031338 0.0008985 -3.488 0.000502 ***
## vowel_ord.L -0.5630529 0.0180089 -31.265 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(time) 6.947e+00 7.273e+00 8.831 5.32e-11 ***
## s(time):frequency 7.789e+00 8.152e+00 121.021 < 2e-16 ***
## s(time):vowel_ordBAT 8.534e+00 8.738e+00 99.586 < 2e-16 ***
## s(frequency) 3.215e-05 3.230e-05 0.241 0.998
## ti(time,frequency) 1.000e+00 1.000e+00 146.535 < 2e-16 ***
## ti(time,frequency):vowel_ordBAT 4.978e+00 5.072e+00 5.354 6.34e-05 ***
## s(time,word_ord) 6.989e+02 9.930e+02 3723.126 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Rank: 1069/1072
## R-sq.(adj) = 1 Deviance explained = 100%
## fREML = -7657.5 Scale est. = 7.464e-06 n = 2100
Three vowels, F1, F2
lexicon_4 <- lexicon_4 %>%
mutate(
vowel_ord = ordered(vowel, levels = c("BART", "BAT", "BET")),
word_ord = as.ordered(word)
) %>%
arrange(word, time) %>%
create_event_start("word")
lexicon_4_gam <- bam(
f1_mean ~
frequency +
vowel_ord +
s(time, bs = "cr") +
s(time, bs = "cr", by = frequency) +
s(time, bs = "cr", by = vowel_ord) +
s(frequency, bs = "cr") +
ti(time, frequency) +
ti(time, frequency, by = vowel_ord) +
s(time, word_ord, bs = "fs"),
data = lexicon_4,
method = "fREML"
)
summary(lexicon_4_gam)
##
## Family: gaussian
## Link function: identity
##
## Formula:
## f1_mean ~ frequency + vowel_ord + s(time, bs = "cr") + s(time,
## bs = "cr", by = frequency) + s(time, bs = "cr", by = vowel_ord) +
## s(frequency, bs = "cr") + ti(time, frequency) + ti(time,
## frequency, by = vowel_ord) + s(time, word_ord, bs = "fs")
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.449284 0.046255 117.811 < 2e-16 ***
## frequency -0.008690 0.002599 -3.343 0.00085 ***
## vowel_ord.L -1.244137 0.022168 -56.123 < 2e-16 ***
## vowel_ord.Q -0.115654 0.022265 -5.194 2.36e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(time) 7.808 7.882 56.477 < 2e-16 ***
## s(time):frequency 4.464 5.050 221.389 < 2e-16 ***
## s(time):vowel_ordBAT 8.740 8.870 142.366 < 2e-16 ***
## s(time):vowel_ordBET 8.699 8.848 141.703 < 2e-16 ***
## s(frequency) 3.832 3.837 4.560 0.00265 **
## ti(time,frequency) 2.936 2.966 112.574 < 2e-16 ***
## ti(time,frequency):vowel_ordBAT 2.534 2.566 4.350 0.01134 *
## ti(time,frequency):vowel_ordBET 2.897 2.929 5.698 0.00257 **
## s(time,word_ord) 669.372 990.000 3862.640 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Rank: 1095/1098
## R-sq.(adj) = 1 Deviance explained = 100%
## fREML = -7756.3 Scale est. = 6.9601e-06 n = 2100
fvisgam(lexicon_4_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BART"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 34.000000.
## * vowel_ord : factor; set to the value(s): BART.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

fvisgam(lexicon_4_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BAT"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 34.000000.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

fvisgam(lexicon_4_gam, view = c("time","frequency"),
cond = list(vowel_ord = "BET"),
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; with 30 values ranging from 1.000000 to 34.000000.
## * vowel_ord : factor; set to the value(s): BET.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##

plot_smooth(
lexicon_4_gam,
view = "time",
cond = list(frequency = 10, vowel_ord = "BAT"),
rug = F,
col = "red",
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; set to the value(s): 10.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##
plot_smooth(
lexicon_4_gam,
view = "time",
cond = list(frequency = 20, vowel_ord = "BAT"),
rug = F,
col = "blue",
add = T,
rm.ranef = TRUE
)
## Summary:
## * frequency : numeric predictor; set to the value(s): 20.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##
plot_smooth(
lexicon_4_gam,
view = "time",
cond = list(frequency = 30, vowel_ord = "BAT"),
rug = F,
col = "green",
add = T,
rm.ranef = TRUE
)

## Summary:
## * frequency : numeric predictor; set to the value(s): 30.
## * vowel_ord : factor; set to the value(s): BAT.
## * time : numeric predictor; with 30 values ranging from 0.000000 to 200000.000000.
## * word_ord : factor; set to the value(s): 1. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(time,word_ord)
##